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Abstract 

This  paper  reviews  the  approach  taken  in  the  literature  to  model  the  effective  transport  coefficients  -  mass  diffusivity,  electrical  conductivity, 
thermal  conductivity  and  hydraulic  permeability  -  of  carbon-fibre  based  porous  electrode  of  polymer  electrolyte  membrane  fuel  cells  (PEMFCs). 
It  is  concluded  that  current  PEMFC  model  do  not  account  for  the  inherent  anisotropic  micro  structure  of  the  fibrous  electrodes.  Simulations 
using  a  2-D  PEMFC  cathode  model  show  that  neglecting  the  anisotropic  nature  and  associated  transport  coefficients  of  the  porous  electrodes 
significantly  influences  both  the  nature  and  the  magnitude  of  the  model  predictions.  This  emphasizes  the  need  to  appropriately  characterize  the 
relevant  anisotropic  properties  of  the  fibrous  electrode. 

©  2006  Elsevier  B.Y.  All  rights  reserved. 
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1.  Introduction 

The  electrode  of  a  polymer  electrolyte  membrane  (PEM)  fuel 
cell  is  a  thin  porous  layer  that  connects  the  flow-field  plate  and 
the  catalyst  layer  as  shown  in  Fig.  1.  This  layer  is  commonly 
referred  to  in  the  literature  as  the  gas  diffusion  layer  (GDL). 
This  term  is  restrictive  and  perhaps  misleading  in  that  it  pre¬ 
supposes  the  primary  function  of  the  porous  electrode  to  be  gas 
diffusion  when  this  is  not  the  case.  Accordingly,  herein  we  use 
the  term  porous  transport  layer  (PTL)  to  describe  the  porous 
electrode  exclusive  of  the  catalyst.  The  most  commonly  em¬ 
ployed  PTL  for  PEM  fuel  cells  are  carbon  cloth  and  carbon 
paper.  These  materials  are  porous  to  allow  for  gas  transport  and 
electronically  conducting  to  ensure  electrons  are  transported  to 
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or  from  the  catalyst  sites.  In  addition,  the  pores  of  the  PTL 
provide  pathway  for  removal  of  liquid  water  generated  in  the 
cathode  catalyst  layer.  Furthermore,  the  PTL  helps  to  conduct 
heat  to  and  from  the  catalyst  layer  and  flow-field  plates.  Clearly, 
the  PTL  serves  as  a  medium  for  transport  of  several  different 
quantities. 

In  recent  years  there  has  been  an  increasing  need  to  engi¬ 
neer  the  PTL  such  that  some  or  all  of  the  transport  processes 
discussed  above  could  be  controlled.  Mathematical  modeling  is 
an  indispensible  tool  for  exploring  new  architectures  and  de¬ 
signs  of  the  PTL.  However,  appropriate  characterization  of  the 
transport  coefficients  of  a  PTL  is  essential  for  meaningful  pre¬ 
dictions  from  such  models.  Over  the  past  20  years,  numerous 
publications  have  been  dedicated  towards  the  modeling  of  poly¬ 
mer  electrolyte  membrane  (PEM)  fuel  cells.  In  a  recent  compre¬ 
hensive  review  of  over  100  papers  on  mathematical  modeling  of 
PEM  and  direct  methanol  fuel  cells  [1],  we  found  that  in  a  ma¬ 
jority  of  previous  work,  the  transport  of  chemical  species  and 
electrons  in  the  PTL  were  treated  simplistically.  For  instance, 
most  of  the  models  for  effective  diffusive  transport  of  chemical 
species  and  effective  electronic  conductivities  in  the  PTL  were 
based  on  correlations  originally  derived  for  granular  porous  me- 
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Fig.  1 .  Section  of  a  PEMFC  cathode  (dotted  line  represents  the  solution  domain 
for  2-D  model). 


dia.  Unfortunately,  the  PTL  micro  structure  does  not  resemble 
such  a  structure.  A  scanning  electron  micrograph  of  the  carbon 
paper  based  PTL  shown  in  Fig.  2  reinforces  this  point.  The  mi¬ 
crostructure  of  the  PTL  can  be  described  as  randomly  laid  out 
cylindrical  fibers  7-10  jam  in  diameter  and  several  millimeters 
long.  Usually  hydrophobic  material  such  as  PTFE  is  added  to 
the  PTL,  apparently  to  avoid  liquid  water  flooding  in  the  elec¬ 
trode  [2].  It  is  clear  that  even  without  the  added  complexity 
of  liquid  water  transport,  the  fibrous  nature  of  PTL  introduces 
anisotropy  with  respect  to  electron  and  gas  transport.  Thus,  the 
commonly  employed  effective  transport  coefficients  originally 
derived  for  granular  porous  media  may  not  be  appropriate  for 
detailed  PEMFC  modeling. 

A  description  of  electron  and  chemical  species  transport  that 
is  representative  of  the  physics  of  the  process  is  essential  because 
the  local  electrochemical  reaction  rates  in  the  catalyst  layers  are 
a  function  of  local  species  concentration  and  the  local  over¬ 
potential  (related  to  electronic  and  ionic  potentials),  which  are 
strongly  influenced  by  transport  through  the  PTL.  In  addition, 
it  is  known  that  convective  transport  does  occur  in  the  PTL  [3], 
and  so  far,  this  directional  dependance  of  the  permeability  has 
not  been  adequately  addressed. 

The  objective  of  this  paper  is  to  review  the  approach  taken  to 
define  the  PTL  effective  transport  coefficients.  Further,  using  a  2- 
Dimensional  model  for  PEMFC  cathode  we  show  that  the  choice 
of  correlations  has  a  significant  impact  on  the  prediction  of  local 
species  fluxes,  species  composition  and  current  densities.  It  is 
important  to  emphasize  that  most  of  the  PEMFC  models  even 
with  an  incorrect  transport  coefficient  can  adequately  fit  selective 
experimental  data  (usually  a  polarization  curve).  However,  such 
models  may  not  be  predictive  and  likely  not  provide  insight  into 
the  influence  of  component  design  and  properties  on  fuel  cell 
performance. 

2.  Generalized  potential  driven  transport 

The  four  key  transport  processes  occurring  in  the  PTL  are 
chemical  species  diffusion,  electron  conduction,  thermal  con¬ 
duction  and  mass  convection.  Each  of  these  processes  are  driven 
by  the  relevant  specific  potential.  The  flux  associated  with  each 
of  these  processes  depend  on  the  gradient  of  the  relevant  poten¬ 


(b) 


Fig.  2.  Typical  carbon  fibre  paper  used  as  porous  transport  layer  in  PEM  fuel 
cells,  (a)  Surface  view  (b)  Edge  view. 


tial  and  effective  transport  coefficient. 

Flux  =  r  effV</>  (1) 

Table  1  summarizes  the  relevant  potential  and  effective  trans¬ 
port  coefficients  for  the  four  processes  occurring  in  the  PTL.  It 
is  useful  to  point  that  the  transport  coefficient  in  itself  depends 
on  the  inherent  property  of  the  transport  medium  and/or  the 
transported  species.  However,  these  transport  coefficients  are 
modified  to  account  for  the  effects  of  the  microstructure  of  the 
transport  media,  which  in  the  present  context  is  the  PTL. 


Table  1 


Key  transport  processes  in  PTL,  associated  driving  potential  and  transport 
coefficient 


Flux 

Transport  coefficient 

Potential 

ni  (molm-2  s-1) 

Am  (m2  s-1) 

Ci  (molm-3) 

q(  Wm“2) 

k  (W  m-1K-1) 

T  (K) 

k  (Am-2) 

(je  (S  m-1) 

0e  (V) 

ip  (Am-2) 

dp  (Sm-1) 

</>P  (V) 

v  (ms-1) 

kill  (m2  Pa-1  s-1) 

P  (Pa) 
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3.  Effective  chemical  species  diffusivity  models 

One  of  the  main  functions  of  the  PTL  is  to  allow  species 
transport  to  and  from  the  catalyst  surface.  Diffusion  of  chemical 
species  is  one  of  the  key  modes  of  transport  of  reactants  from  the 
flow-channels  to  the  catalyst  layer.  For  diffusion  through  porous 
material,  the  binary  diffusion  coefficient,  D^,  is  modified  to 
account  for  the  effects  of  the  microstructure  of  the  porous  media: 

Df  =  /(e)  •  Dij  (2) 

It  must  be  noted  that  all  microstructural  effects  are  lumped 
in  a  function  f(e)  of  a  single  parameter — porosity  (e).  Various 
theories  are  available  to  determine  the  function  f(e). 

The  Bruggeman  relation  is  an  analytic  expression  which  be¬ 
longs  to  the  class  of  solutions  referred  to  as  Effective  Medium 
Approximations.  The  chemical  species  diffusivity  in  the  PTL  is 
most  commonly  described  by  the  Bruggeman’s  relation,  wherein 
the  characteristic  function  describing  the  porous  micro  structure 
is  given  by: 

m  =  e1'5  (3) 


Fig.  3.  Comparison  of  effective  diffusivity  from  percolation  theory  with  the 
Bruggeman  correlation. 


Interestingly,  this  correlation  was  originally  derived  for  effec¬ 
tive  dispersion  medium  permittivity  of  beds  of  spherical  particles 
of  differing  sizes  [4]  and  was  later  verified  experimentally  by 
measuring  the  conductivity  of  an  electrolyte  surrounding  spher¬ 
ical  glass  beads.  The  relation  fits  experimental  data  for  diffusiv¬ 
ity  in  granular  packed  beds  [5].  However,  it  appears  that  when 
this  correlation  was  first  adopted  to  describe  mass  diffusivity  in 
fibrous  PTL,  the  validity  of  the  correlation  applied  to  fibrous 
porous  structures  was  not  considered. 

3.1.  Diffusivity  of  gases  in  fibrous  porous  media 


Except  for  a  paper  by  Nam  and  Kaviany  [6],  none  of  the 
PEMFC  models  have  considered  the  anisotropy  of  fibrous  PTLs 
with  respect  to  species  diffusivity.  It  is  interesting  to  note  that  a 
considerable  amount  of  work,  both  theoretical  and  experimen¬ 
tal,  has  been  undertaken  to  understand  and  describe  the  phenom¬ 
ena  of  gas  diffusivity  in  fibrous  porous  media.  Tomadakis  and 
Sotirchos  [7,8]  applied  Monte  Carlo  simulation  to  estimate  the 
effective  diffusivity  of  chemical  species  in  a  packing  of  overlap¬ 
ping  fibres  randomly  distributed  in  a  plane.  Nam  and  Kaviany 
[6]  reported  the  correlation  given  by  Eq.  (4),  originally  reported 
by  Tomadakis  and  Sotirchos  [8],  where  a  =  0.521  for  the  in¬ 
plane  direction,  a  =  0.785  for  the  through-plane  direction,  and 
£p  =  0.11. 


f(s)  =  £ 


Fig.  3  compares  the  Bruggeman  relation  to  the  anisotropic 
percolation  type  correlation  of  Tomadakis  and  Sotirchos.  It  is 
to  note,  that  with  this  model,  the  through-plane  value  is  always 
lower  than  the  in-plane  value,  and  that  the  differences  decrease 
with  increasing  porosity.  At  a  porosity  of  0.8,  the  in-plane  value 
is  only  91%  of  the  value  suggested  by  the  Bruggeman  rela¬ 
tion,  while  at  a  porosity  of  0.4,  the  in-plane  value  is  only  65% 
of  the  Bruggeman  relation.  This  indicates  that  when  diffusion 


is  an  important  mechanism  the  Bruggeman  relation  will  over¬ 
predict  the  concentration  of  Oxygen  in  the  catalyst  layer  and 
hence  the  current  density.  The  extent  of  this  effect  is  analyzed  in 
Section  8. 


4.  Thermal  conductivity 


Contrary  to  some  of  the  modeling  results  [9],  experimental 
data  [10]  indicate  that  significant  temperature  gradients  exist 
across  the  PTL.  One  likely  reason  that  the  models  are  unable  to 
capture  experimentally  observed  temperature  gradients  is  that 
inaccurate  and  unrealistically  large  values  for  thermal  conduc¬ 
tivity  of  the  PTL  is  employed  in  the  model.  In  our  survey  of  over 
100  fuel  cell  modelling  papers,  every  non-isothermal  model  used 
the  volume  averaged  thermal  conductivity.  A  physical  interpre¬ 
tation  of  this  volume-averaging  is  that  of  a  porous  media  that  is 
perfectly  sorted  such  that  all  the  solid  material  is  parallel  to  all 
the  void  space,  as  shown  in  Fig.  4.  This  means  that  some  of  the 
heat  energy  is  forced  either  to  pass  through  the  solid  phase,  and 
some  is  forced  through  the  void  phase  but  no  path  involving  a 
combination  of  solid  and  void  is  allowed. 


^eff,  parallel  —  +  (1  o)ks 


Another  possible  way  of  sorting  the  material  is  to  stack  re¬ 
gions  of  solid  and  regions  of  void  such  that  all  of  the  heat  energy 
is  forced  to  pass  through  both  solid  and  void.  This  ordering  cor¬ 
responds  to  a  series  resistance  network,  and  results  in  a  volume 
weighted  harmonic  average  for  the  effective  conductivity. 


, series  — 


ksk y 

sks  +  (1  —  s)ky 


These  two  expressions  are  plotted  in  Fig.  4  for  the  case  of 
ks/kv  =  100,  along  with  the  series  and  parallel  resistance  net¬ 
works  from  which  they  are  derived.  The  parallel  resistance  pro¬ 
vides  an  upper  limit  for  the  effective  conductivity,  while  the 
series  resistance  network  provides  a  lower  limit.  Clearly,  real 
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Fig.  4.  Normalized  thermal  conductivity  for  series  ordering  and  parallel  order¬ 
ing.  ks/kw  =  100. 

porous  media  with  random  orderings  of  solid  and  void  space 
must  lie  somewhere  in  between  the  two  limits. 

In  a  recently  completed  thesis  [11]  the  anisotropic  thermal 
conductivity  of  carbon  fibre  mattes  was  computed  using  a  dis¬ 
cretization  of  fibrous  microstructures.  The  results  indicate  that 
the  volume  average,  or  parallel  resistance  model,  significantly 
over-predicts  the  actual  thermal  conductivity  which  lies  closer  to 
the  series  approximation — especially  in  the  in-plane  direction. 
This  means  that  predicted  temperature  gradients  across  the  PTL 
will  be  significantly  lower  than  in  reality. 

5.  Electronic  conductivity 

Although  several  models,  citing  high  electronic  conductivity 
of  the  medium,  ignore  the  electronic  potential  gradients  in  the 
PTL  it  has  been  shown  [12-14]  that  the  electron  transport  in  the 
PTL  can  influence  the  distribution  of  current  density  or  reaction 
rate  in  the  catalyst  layer.  The  effective  electronic  conductivity  of 
PTL  materials,  to  a  first  approximation,  is  analogous  to  the  ef¬ 
fective  thermal  conductivity  in  the  limit  of  the  void  conductivity 
taken  to  zero.  For  series  arrangement,  the  effective  conductivity 
is  zero  because  the  inter-layer  of  void  serves  as  a  perfect  insu¬ 
lator.  The  parallel  arrangement  yields  the  following  expression 
for  effective  conductivity: 

^eff,  parallel  =  (1  —  o)(7s  (7) 

These  approximations  are  seldom  used  and  effective  conductiv¬ 
ity  values  instead  are  referenced  from  other  papers  or  manufac¬ 
turers  web  sites.  It  is  only  recently  that  the  effect  of  anisotropic 
electronic  conductivity  has  been  investigated  [14,15].  While 
there  is  very  little  experimental  verification  of  anisotropic  elec¬ 
tronic  conductivities,  some  evidence  seems  to  suggest  that  the 
in-plane  conductivity  is  about  an  order  of  magnitude  higher  than 
the  through  plane  value  [2].  It  will  be  shown  in  section  8  that 
higher  in-plane  conductivities  can  drastically  alter  the  location 
of  current  production  by  electrochemical  reaction. 


Fig.  5.  The  effect  of  anisotropic  permeability  on  channel  pressure  drop.  Re  = 
100,  L  =  40  mm. 

6.  Hydraulic  permeability 

The  permeability  of  PTL  material  is  of  most  significance 
when  serpentine  flow  fields  are  used  because  this  induces  an 
in-plane  pressure  gradient  on  the  porous  transport  media.  This 
pressure  gradient  results  in  the  diversion  of  significant  flow  un¬ 
der  the  lands  and  causes  a  substantial  decrease  in  the  pressure 
drop  required  to  drive  the  flow.  This  ‘feature’  has  been  demon¬ 
strated  numerically  [3],  and  has  been  exploited  physically  to 
detect  flooding  in  the  porous  transport  layer  [16]. 

It  can  be  expected  that  similar  to  other  properties,  the  perme¬ 
ability  of  PTL  materials  is  also  likely  to  be  anisotropic,  though 
these  values  are  not  often  reported.  The  relative  importance  of  in¬ 
plane  and  through-plane  permabilities  was  recently  considered 
[3],  and  it  was  shown  that  the  in-plane  value  is  most  significant. 
This  is  demonstrated  in  Fig.  5  which  plots  the  pressure  drop 
necessary  to  drive  a  fixed  flow  rate  through  a  single  serpentine 
channel  mounted  on  PTLs  of  varying  permeability.  In  the  case  of 
isotropic  PTLs  the  pressure  drop  follows  a  sigmoidal  type  curve 
with  a  low  pressure  drop  corresponding  to  high  permeability  and 
a  much  larger  pressure  drop  in  the  limit  of  impermeable  PTLs. 
In  this  case,  an  isotropic  PTL  is  virtually  impermeable  when  the 
permeability  is  of  order  10-12  m2.  Interestingly,  if  the  through- 
plane  permeability  is  held  constant  while  the  in-plane  perme¬ 
ability  is  increased,  the  driving  pressure  gradient  decreases  by 
over  40%.  Unfortunately,  virtually  all  measurements  of  PTL  per¬ 
meabilities  correspond  to  the  through-plane  values.  This  can  be 
very  misleading,  and  accordingly,  there  is  a  significant  need  to 
characterize  PTLs  for  in-plane  permeability. 

7.  Porosity  of  the  porous  transport  layer 

As  pointed  out  earlier,  correlations  describing  effective  trans¬ 
port  coefficients  account  for  porous  media  effects  in  terms  of 
a  single  parameter — porosity.  As  such,  porosity  is  an  impor- 
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Fig.  6.  Picture  of  a  PTL  showing  localized  compression  under  the  land  area  [17] 
(dark  striations  correspond  to  PTL  under  the  land  area  that  have  been  subjected 
to  compression). 

tant  model  input  parameter.  Selecting  a  representative  value  of 
PTL  porosity  is  complicated  because  in  an  operational  fuel  cell 
the  PTL  is  subjected  to  spatially  variant  stress  or  compression, 
which  in  turn  results  in  a  spatially  distributed  porosity.  For  in¬ 
stance,  the  PTL  under  the  land  area  is  expected  to  be  compressed 
to  a  greater  extent  than  that  under  the  channel,  consequently  PTL 
porosity  is  expected  to  be  lower  than  that  under  the  channel  area. 
This  effect  is  highlighted  in  Fig.  6,  which  shows  a  picture  of  a 
carbon-fibre  based  PTL  that  had  been  placed  between  two  plates, 
each  12.7  cm  x  12.7  cm  -  one  which  was  flat  and  another  plate 
with  single  serpentine  channel  -  subjected  to  a  stress  of  1 .2  MPa 

[17] .  The  dark  striations  correspond  to  the  PTL  that  were  under 
the  land  are  were  subjected  to  compression.  The  extent  to  which 
PTL  porosity  is  compressed  under  various  stress  or  mechanical 
load  conditions  has  been  reported  by  [18].  The  extent  of  com¬ 
pression  or  more  correctly  strain  was  found  to  vary  non-linearly 
with  the  applied  stress.  The  porosity  of  uncompressed  PTL  can 
be  over  80%  [19,18]  but  as  per  data  of  Mishra  and  co-workers 

[18]  it  can  reduce  to  55-65%  when  subjected  to  local  stress  of 
0.5-1  MPa.  In  contrast,  much  of  the  original  PEMFC  modeling 
work  have  employed  a  porosity  value  of  40%.  In  the  follow¬ 
ing  section,  we  examine  the  extent  to  which  the  PTL  porosity 
influences  performance  of  a  PEMFC  via  simulations. 

8.  Effect  of  transport  coefficients  on  the  predictions  of 
fuel  cell  performance 

So  far,  we  have  discussed  the  applicability  of  various  trans¬ 
port  coefficients  and  the  need  for  anisotropic  transport  coeffi¬ 
cients  on  physical  grounds.  We  will  now  demonstrate  the  im¬ 
pact  of  transport  coefficients  and  their  directional  dependance 
on  fuel  cell  performance  using  a  simplified  two-dimensional  (2- 
D)  PEMFC  cathode  model  applied  in  the  cross-channel  plane. 


chemical  processes:  multi-component  diffusion  of  chemical 
species  in  the  PTL  and  catalyst  layer,  oxygen  dissolution  and 
diffusion  within  idealized  catalyst  agglomerates,  electron  con¬ 
duction  in  the  PTL  and  catalyst  layer,  proton  transport  in  the 
catalyst  layer  and  electrochemical  reaction  governed  by  Butler- 
Volmer  type  kinetics. 

The  problem  domain,  which  assumes  symmetry  at  the  chan¬ 
nel  centreline  and  the  land  centreline,  is  shown  in  Fig.  1 .  It  is  use¬ 
ful  to  note  that  this  symmetry  is  not  strictly  valid  for  the  case  of 
convective  transport  under  the  lands,  a  situation  not  investigated 
in  our  2-D  model.  The  coupled  transport-reaction  equations  are 
solved. 


8.1.1.  Species  transport 

The  transport  of  chemical  species  in  the  PTL  and  the  catalyst 
layer  is  modeled  as  multi-component  diffusion  in  the  porous 
media  by  the  Maxwell-Stefan  equation: 


V  •  Hi  =  V  • 


(vwj  +  Wj 


where,  ni  is  mass  flux  (kgm~2  s-1),  p  gas  mixture  density,  and 
Wi  is  mass  fraction  of  species  i.  The  average  molecular  weight 
of  gas  mixture  (M)  is: 


M  =  'ffxiMi 


and  the  gas  mixture  density  is  is  calculated  by: 


P tot 

RTo 


(10) 


In  the  catalyst  layer,  the  gradient  of  oxygen  molar  flux  (No2) 
is  equal  to  the  rate  of  reaction: 


V  •  Nq2  =  —  Rq2 


(11) 


The  rate  of  reaction,  Ro2 ,  within  the  catalyst  layer  is  described 
by  the  following  equation: 


1  +  (r  agg  +  <$)<$ 

Erkc(l  —  £cat)  ^agg^agg^  _ 


(12) 


The  above  Eq.  (12)  describes  the  rate  of  oxygen  consumption 
in  terms  of  the  oxygen  concentration  at  the  agglomerate  surface 
(Co9  J,  an  reaction  effectiveness  factor,  Ev ,  to  account  for  the 
diffusion-reaction  process  within  the  agglomerate,  and  an  elec¬ 
trochemical  rate  constant,  kc.  The  reaction  effectiveness  factor 
for  spherical  agglomerate  is  calculated  by: 


Er  =  T  ( — - L) 

<Pl  \tanh(30L)  30l  / 

where 


(13) 


8.1.  2D  PEMFC  cathode  model 

The  details  of  the  2-D  model  formulation  are  given  in  [  1 3 , 1 5] . 
The  model  accounts  for  the  following  key  physico-electro- 


ragg 

3 


kc 

£>eff 


Dc  ff  =DsL5 


(14) 

(15) 
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The  electrochemical  reaction  rate,  kc ,  is  essentially  Butler- 
Volmer  equation  accounting  for  catalyst  micro  structure  effects: 


(  ^Pt^ac  \ 

1 

4-h 

CD 

U 

O 

•  ^ 

i _ 

\  4Ffc](l  —  £c  at)/ 

-  ^C>2, ref  - 

x  [exp (acfri)  -  exp((l  -  ac)fy)] 

(16) 

where,  /  =  F/RT  and  Sac  =  3/(rptppt).  All  other  parameters  in 
the  above  Eq.  16  are  reported  in  Table  3  and  the  variable  77,  the 
overpotential  is  an  input  variable.  More  details  of  the  catalyst 
layer  model  are  provided  in  [13]. 

8.1.2.  Charge  transport 

The  electron  transport  within  the  PTL  is 
law: 

described  by  Ohm’s 

4,  PTL  =  -kfpTL^Cps 

(17) 

also, 

V  •  /e,  PTL  =  0 

(18) 

The  effective  conductivity  of  the  PTL  is  specified  as  an  input 
parameters  for  isotropic  and  anisotropic  models. 

The  electron  and  proton  transport  within  the  catalyst  layer 
are  described  by  the  following  equations: 

4,  cat  —  ^sfcat^^s 

(19) 

ip,  cat  =  —  ^l^at^^l 

(20) 

also, 

^ (4,  cat  T  ip,  cat)  =  0 

(21) 

and 

^  ’  4,  cat  —  V  •  i 

(22) 

The  effective  conductivities  of  the  proton  and  electron  are 
calculated  by  Bruggeman  relation  because  the  catalysts  are  con¬ 
sidered  to  be  spherical  agglomerates: 


Table  2 


Operating  parameters  considered  in  2-D  model 


Inlet  pressure  (absolute)  (atm) 

1.5 

Operating  temperature  (°C) 

80 

Inlet  oxygen/nitrogen  ratio 

21/79 

Inlet  air  relative  humidity 

50 

Table  3 

Input  geometric  and  microstructural  parameters 

Parameters  (units) 

Value 

PTL  thickness  ([xm) 

300 

Porosity  of  PTL,  £ptl 

0.6 

Pt  loading,  mpt  (mg  cm-2  Pt) 

0.4 

Pt  particle  diameter,  rpt  (nm) 

2.5 

Radius  of  agglomerate,  ragg  (|xm) 

1 

Effective  specific  agglomerate  surface  area,  aa gg  (m2  m-3) 

3.6  X  105 

Catalyst  layer  thickness,  tc\  (fxm) 

15 

Porosity  of  catalyst  layer,  scat 

0.1 

Electrolyte  film  covering  each  agglomerate,  8  (nm) 

80 

8.1.3.  Boundary  conditions 

For  species  transport,  no  flux  boundary  conditions  was  spec¬ 
ified  at  the  right,  left  and  bottom  sides  of  the  domain  in 
Fig.  1.  At  the  top  side  of  the  domain,  concentrations  of  oxygen, 
nitrogen  and  water  vapor  was  specified  at  the  interface  between 
the  PTL  and  the  channel,  whereas  no  flux  boundary  condition 
was  implemented  at  the  the  interface  between  PTL  and  the  land 
region.  Similarly,  for  electron  transport  no  flux  boundary  con¬ 
ditions  was  specified  at  the  right,  left  and  bottom  sides  of  the 
domain.  At  the  top  side  of  domain,  no  flux  boundary  condition 
was  implemented  at  the  the  interface  between  the  PTL  and  the 
channel,  whereas  an  electronic  phase  potential  of  zero  was  spec¬ 
ified  at  the  interface  between  the  PTL  and  the  land.  For  proton 
transport,  again,  no  flux  boundary  conditions  were  specified  at 
the  right  and  left  hand  sides  of  the  domain.  The  bottom  side  of 
the  domain  was  set  to  proton  phase  potential  equal  to  the  spec¬ 
ified  overpotential  whereas  the  interface  between  the  PTL  and 
the  catalyst  layer  was  set  as  no  flux  boundary. 

The  key  model  input  parameters  are  provided  in  Tables  2-4. 
The  details  of  how  these  parameters  were  obtained  are  provided 


fcfcat  =  JfcsKl  -  ScatXl  -  fiagg)]15  (23) 

and 

^Uat  =  Md  -  ScatKgg]1'5  (24) 


The  proton  conductivity  of  the  electrolyte  membrane  is  rep¬ 
resented  as  that  of  the  Nation  and  modeled  in  [13]  according  to 
the  following  equations  based  on  work  by  Springer  group  [20] : 


k\  =  100[0.005139A 


0.00326]  exp 


1268 


and,  the  local  water  content,  A,  in  the  electrolyte  phase  depends 
on  the  relative  humidity  (RH): 

A  =  0.3  +  10.8RH- 16RH2  +  14.1RH3  (26) 


Table  4 


Physical  and  electrochemical  kinetic  parameters 


Parameters  (units) 

Value 

Conductivity  of  PTL,  ks  (S  m-1) 

100 

O2  Diffusivity  in  Nafiona,  D  (m2s_1) 

8.45  x  nr10 

Binary  diffusivities,  Dij  (atm  m2s_1 ) 

Ad2-h2o 

3.70  x  10“5 

£>o2-n2 

2.79  x  10“5 

^n2-h2o 

3.87  x  10“5 

Cathodic  transfer  coefficient,  ac 

Low  slope  (>0.8  V) 

1 

High  slope  (  <0.8  V) 

0.61 

Reference  exchange  current  densitya,  z'o  (A  cm-2) 

Low  slope  (>0.8  V) 

3.85  x  10-8 

High  slope  (<0.8  V) 

1.5  x  10“6 

Reference  O2  concentration13,  Cq2  ref(molm~3) 

0.85 

a  80  °C. 

b  80  °C  and  1.5  atm. 
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in  [13].  The  governing  equations  with  their  associated  boundary 
conditions  were  solved  using  the  Multiphysics  software  FEM- 
LAB  (Comsol).  The  simulation  results  are  presented  in  the  fol¬ 
lowing  sub-sections.  The  effects  of  anisotropic  mass  diffusivity 
and  electronic  conductivity  on  the  polarization  curve  and,  more 
importantly,  on  the  current  density  distribution  in  cross-channel 
direction  are  investigated. 

8.2.  Isotropic  and  anisotropic  mass  diffusivity 

In  this  section  we  examine  the  influence  of  mass  diffusiv¬ 
ity  magnitude  by  comparing  three  different  effective  diffusivity 
models  at  two  different  PTL  porosity  values.  The  three  mod¬ 
els  are:  (i)  isotropic  effective  diffusivity  calculated  by  Brugge- 
man  relation,  (ii)  isotropic  effective  diffusivity  calculated  using 
in-plane  percolation  model,  and  (iii)  anisotropic  diffusivity  cal¬ 
culated  by  percolation  model.  The  anisotropic  model  consid¬ 


ers  variation  of  diffusivity  in  two  directions  only — in-plane  and 
through-plane. 

8.2.1.  Mass  diffusivity  effects  for  low  porosity  PTL 

We  first  present  the  simulation  results  for  the  case  of  PTL  with 
a  low  porosity  of  0.40,  which  is  likely  an  unrealistic  value  but  one 
that  has  been  used  commonly  in  PEMFC  modeling  literature. 
Fig.  7  presents  simulation  results  for  three  different  effective 
diffusivity  models.  Fig.  7(a  and  b)  depicts  the  distribution  of 
current  at  the  membrane  catalyst  interface  for  two  different  load 
conditions,  A  and  B  indicated  on  cathode  polarization  curve, 
Fig.  7(c). 

If  we  compare  the  two  isotropic  models,  it  can  be  observed 
that  at  both  load  conditions,  the  predicted  current  is  everywhere 
lower  in  the  case  of  the  percolation  model  because  of  the  lower 
effective  mass  diffusivity.  This  effect  is  more  pronounced  under 
the  land  area.  The  maximum  decrease  is  on  the  order  of  8%  at 


0.2  0.4  0.6  0.8  1 

(c)  Current  Density  [A/cm2] 


Fig.  7.  Effect  of  diffusivity  models  on  predicted  (a-b)  current  density  distributions  and  (c)  polarization  curves  for  PTL  porosity  of  0.4 
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the  lower  load  condition  and  16%  at  the  higher  load  condition. 
It  is  also  evident  from  the  polarization  curve,  Fig.  7(c)  that  the 
onset  of  mass  transfer  limitations  occurs  at  significantly  lower 
currents  with  decreasing  effective  diffusivity.  It  is  of  note  that 
the  model  employed  to  investigate  these  effects  is  a  single  phase 
model  such  that  there  is  no  liquid  water  present  in  any  of  the 
simulations.  The  sharp  voltage  drop  off  at  higher  currents  can 
be  attributed  to  a  mass  transfer  resistance  through  the  Nation 
film  which  is  assumed  to  surround  each  catalyst  agglomerate 
[13].  The  thickness  of  the  Nation  film  has  been  held  constant  for 
all  cases  considered  herein. 

Next,  we  consider  the  effect  of  anisotropic  mass  diffusiv¬ 
ity,  with  both  the  in-plane  and  through-plane  diffusivities  deter¬ 
mined  from  the  percolation  theory  (Eq.  (4).  If  we  compare  the 
isotropic  percolation  and  anisotropic  percolation  models,  it  can 
be  observed  from  Fig.  7(a  and  b)  that  there  is  a  significant  in¬ 
crease  in  the  predicted  current  under  the  land,  while  the  current 
under  the  channel  remains  virtually  unchanged.  This  is  due  to 


the  fact  that  the  anisotropic  model  results  in  a  relative  increase 
in  the  through  plane  conductivity.  Also  evident  is  a  slight  in¬ 
crease  in  the  maximum  current  predicted,  consistent  with  the 
above  arguments.  An  important  point  to  be  emphasized  is  that 
the  PEMFC  performance  predicted  using  Bruggeman  model  is 
significantly  higher  than  than  anisotropic  model. 

8.2.2.  Mass  diffusivity  effects  for  high  porosity  PTL 

In  this  section  we  examine  to  what  extent  do  the  predictions 
from  the  three  diffusivity  models  vary  for  an  input  PTL  porosity 
of  0.6,  a  value  more  representative  of  a  PTL  in  an  operating 
PEMFC.  Fig.  8 (a  and  b)  present  the  current  density  distribution 
for  two  different  load  conditions,  B  and  C,  indicated  on  cathode 
polarization  curve,  Fig.  8(c). 

In  comparison  to  results  presented  in  preceding  section 
(Fig.  7),  several  points  may  be  noted.  For  PTL  porosity  of  0.6, 
each  of  the  three  models  exhibit  a  limiting  current  density  that  is 
at  least  20%  higher  than  that  for  PTL  porosity  of  0.4  Examining 


Fig.  8.  Effect  of  diffusivity  models  on  predicted  (a-b)  current  density  distributions  and  (c)  polarization  curves  for  PTL  porosity  of  0.6. 
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Figs.  7(c)  and  8(c),  it  can  be  observed  that  differences  in  predic¬ 
tions  from  the  three  models  begin  to  appear  at  Cathode  potential 
of  less  than  0.8  V  and  less  than  0.6  V  for  PTL  porosities  of  0.4 
and  0.6,  respectively.  Comparing  the  current  density  distribu¬ 
tions  in  Figs.  7(b)  and  8(a)  which  correspond  to  load  Point  B,  it 
can  be  concluded  that  upon  increasing  the  PTL  porosity  from  0.4 
to  0.6,  the  facility  of  oxygen  transport  vis-a-vis  electron  trans¬ 
port  is  enhanced.  This  effect  is  manifested  in  a  higher  current 
density  generated  being  under  the  land  than  that  under  channel 
for  PTL  porosity  of  0.6  For  the  case  of  PTL  porosity  of  0.6,  the 
effect  of  mass  transfer  limitations  become  prominent  at  higher 
overpotential  or  lower  cathode  potential  as  observed  in  higher 
current  densities  under  the  land  area  in  Fig.  8(b).  An  underlying 
point  to  be  noted  that  the  commonly  employed  Bruggeman  rela¬ 
tion  predicts  a  higher  current  density  compared  to  that  by  other 
models,  although  the  differences  are  within  5-6%. 

From  the  results  presented  in  Sections  8.2.1  and  8.2.2,  it  is 
clear  that  predicted  fuel  cell  performance  is  sensitive  to  the  ab¬ 


solute  value  of  the  effective  diffusivity  used  in  the  porous  media. 
This  value  is  currently  not  well  represented  by  the  Bruggeman 
relation  originally  developed  for  isotropic  porous  media  com¬ 
posed  of  dense  packings  of  solid  spheres.  It  would  be  a  signifi¬ 
cant  advantage  to  fuel  cell  modelling  if  the  determination  of  the 
effective  diffusivity  was  a  standard  test  in  the  characterization 
of  fuel  cell  PTLs. 

The  local  reaction  rate  is  coupled  to  the  local  electronic  po¬ 
tential,  ionic  potential,  concentration  and  temperature  through 
the  Butler- Volmer  equation.  So  far  we  have  explored  only  the 
effect  of  concentration  changes  resulting  from  different  effec¬ 
tive  diffusivities.  We  will  now  investigate  the  effect  of  electronic 
conductivity. 

8.3.  Isotropic  and  anisotropic  electronic  conductivity 

As  discussed  earlier,  it  is  expected  that  the  electronic  conduc¬ 
tivity  of  fibrous  PTL  materials  will  be  anisotropic  although  most 


Fig.  9.  Effect  of  anisotropic  electronic  conductivity  on  predicted  (a-b)  current  density  distribution  and  (c)  polarization  curves. 
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models  assume  isotropy.  In  order  to  quantify  the  effect  of  this  as¬ 
sumptions  two  cases  are  compared:  one  assuming  isotropic  elec¬ 
tronic  conductivity,  and  the  other  with  an  in-plane  conductivity 
one  order  of  magnitude  higher  than  the  through-plane  value.  In 
order  to  isolate  the  effect  of  changes  in  electronic  conductivity, 
the  mass  diffusivity  is  isotropic  and  described  by  the  Bruggeman 
relation  for  the  two  cases. 

Fig.  9(a  and  b)  presents  the  current  density  distribution  at 
the  catalyst  layer-PTL  interface  for  two  different  load  points  A 
and  B  indicated  on  the  polarization  curve  in  Fig.  9(c).  The  same 
two  load  points  for  the  cases  of  isotropic  electronic  conductivity 
and  anisotropic  electronic  conductivity  was  maintained.  In  the 
case  of  isotropic  electronic  conductivity,  for  load  point  A  the 
reaction  is  limited  by  the  relatively  low  value  of  the  conductivity 
(100  S  m-1)  such  that  the  current  is  higher  under  the  land  area, 
with  the  maximum  occurring  under  the  land  area.  For  load  point 
A,  as  the  current  density  is  increased  the  location  of  maximum 
current  is  increased  to  region  near  the  channel-land  interface. 
For  the  case  of  anisotropic  conductivity,  because  the  in-plane 
electronic  conductivity  is  10-folds  through-plane  conductivity, 
the  resistance  to  electron  transport  from  land  area  to  under  the 
channel  area  is  significantly  reduced,  and  the  maximum  shifts 
to  the  channel  centreline  where  the  oxygen  concentrations  are 
highest. 

This  is  an  extremely  interesting  case,  because  the  predicted 
polarization  curves  are  almost  identical,  and  yet  the  local  current 
distribution  is  completely  different.  Unfortunately,  the  current 
method  for  validation  of  PEMFC  models  are  based  on  compar¬ 
ing  model  predictions  to  experimental  data  for  a  single  polariza¬ 
tion  curve.  This  remarkable  coincidence  illustrates  the  danger 
of  attempting  to  validate  fuel  cell  models  using  a  single  polar¬ 
ization  curve.  Clearly,  the  choice  of  physics  modeled  and  the 
setting  of  the  parameters  in  those  models  has  a  large  impact  on 
the  detailed  results.  These  differences,  while  significant  for  the 
design  and  optimization  of  fuel  cells  are  completely  obscured 
upon  integration  to  the  polarization  curve. 

Fuel  cell  models  require  many  transport  and  reaction  param¬ 
eters,  many  of  which  are  currently  not  well  characterized.  It  is 
crucial  to  recognize  the  physical  basis  for  these  parameters,  and 
to  account  for  the  microstructure  of  the  electrodes  used  in  the 
fuel  cell  being  modeled.  In  some  cases  certain  physics  are  ne¬ 
glected  in  a  model  and  the  one  or  more  remaining  parameters 
are  tuned  to  match  a  single  polarization  curve.  This  will  clearly 
detract  from  the  predictive  capability  of  the  derived  model.  It  is 
our  opinion  that  model  validation  must  involve,  at  a  minimum, 
agreement  with  multiple  polarization  curves  taken  under  differ¬ 
ent  operating  conditions  or  with  different  fuel  cell  components 
relevant  to  the  study  in  question. 

9.  Conclusions 

In  this  paper,  the  approach  taken  in  the  literature  to  describe 
effective  transport  coefficients  of  porous  PEMFC  electrodes 
was  reviewed.  It  was  concluded  that  a  majority  of  the  exist¬ 
ing  PEMFC  models  volume- average  the  effect  of  porous  media 
using  a  function  that  assumes  the  electrode  to  be  isotropic.  This 
is  despite  the  fact  that  the  inherently  anisotropic  microstructure 


of  the  porous  carbon-fibre  electrode  results  in  a  distinctly  differ¬ 
ent  effective  mass  diffusivity,  electronic  conductivity,  thermal 
conductivity,  and  hydraulic  permeability  in  the  through-plane 
and  the  in-plane  directions. 

Using  a  2-D  PEMFC  cathode  model,  it  was  shown  that  treat¬ 
ing  the  electrodes  as  an  isotropic  porous  media  yields  signif¬ 
icantly  different  current  density  predictions  than  anisotropic 
treatments.  Most  importantly,  it  was  demonstrated  that  for  cer¬ 
tain  sets  of  parameters,  both  isotropic  and  anisotropic  models 
yield  virtually  identical  polarization  curves,  however,  the  cur¬ 
rent  density  distributions  are  completely  different.  This  result 
highlights  the  danger  of  validating  detailed  PEMFC  models  with 
a  single  set  of  data,  namely,  the  polarization  curve.  More  impor¬ 
tantly,  it  clearly  identifies  the  need  for  appropriately  character¬ 
izing  the  effective  transport  coefficients  of  the  porous  such  that 
correct  representation  of  the  microstructural  effects  and  physics 
of  the  processes  are  accounted  for.  At  a  minimum,  correct  poros¬ 
ity  of  the  PTL  must  be  determined  under  realistic  conditions.  The 
spatial  variation  of  porosity  must  be  also  of  accounted  for  in  the 
PEM  models. 
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